data5b <- as.data.frame(x5b$summary)
data5b <- data5b %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data5b <- data5b %>%
mutate(time = c(0, 1, 2, 3, 4))
################################################################
# Plots - VSP
################################################################
pdf(file = "Figure2-D.pdf",   # The directory you want to save the file in
width = 9, # The width of the plot in inches
height = 3)
par(mar = c(1.5, 2, 2 , 1), oma = c(4, 4,1.5, 0))
par(mfrow = c(1, 5))
plot(x = data1b$time,
y = data1b$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-60, 60, 6),
xlab = "", ylab = "",
pch = 16, cex.axis=1.25)
segments(x0 = data1b$time, x1 = data1b$time,
y0 = data1b$lower_b, y1 = data1b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data2b$time,
y = data2b$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data2b$time, x1 = data2b$time,
y0 = data2b$lower_b, y1 = data2b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data3b$time,
y = data3b$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data3b$time, x1 = data3b$time,
y0 = data3b$lower_b, y1 = data3b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data4b$time,
y = data4b$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data4b$time, x1 = data4b$time,
y0 = data4b$lower_b, y1 = data4b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data5b$time,
y = data5b$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data5b$time, x1 = data5b$time,
y0 = data5b$lower_b, y1 = data5b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
mtext("Mahalanobis", line = -0.2, at = 0.11, outer = TRUE, cex = 1)
mtext("PS \n Matching", line = -1.8, at = 0.31, outer = TRUE, cex = 1)
mtext("PS \n Weighting", line = -1.8, at = 0.51, outer = TRUE, cex = 1)
mtext("CBPS \n Matching", line = -1.8, at = 0.71, outer = TRUE, cex = 1)
mtext("CBPS \n Weighting", line = -1.8, at = 0.91, outer = TRUE, cex = 1)
mtext("Estimated Effect of VSP", side = 2, line = .3, at = 0.5,
outer = TRUE, cex = 1)
mtext(1, text = "Years relative to the administration of treatment", line = 1.2,
at = 0.5, outer = TRUE, cex = 1)
dev.off()
################################################################
# Preparing the data - Sanctuary
################################################################
x1c <- summary(pe.results5.s2_maha5b)
x2c <- summary(pe.results5.s2_ps.m5)
x3c <- summary(pe.results5.s2_ps.w5)
x4c <- summary(pe.results5.s2_CBPSm5)
x5c <- summary(pe.results5.s2_CBPSw5)
data1c <- as.data.frame(x1c$summary)
data1c <- data1c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data1c <- data1c %>%
mutate(time = c(0, 1, 2, 3, 4))
data2c <- as.data.frame(x2c$summary)
data2c <- data2c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data2c <- data2c %>%
mutate(time = c(0, 1, 2, 3, 4))
data3c <- as.data.frame(x3c$summary)
data3c <- data3c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data3c <- data3c %>%
mutate(time = c(0, 1, 2, 3, 4))
data4c <- as.data.frame(x4c$summary)
data4c <- data4c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data4c <- data4c %>%
mutate(time = c(0, 1, 2, 3, 4))
data5c <- as.data.frame(x5c$summary)
data5c <- data5c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data5c <- data5c %>%
mutate(time = c(0, 1, 2, 3, 4))
################################################################
# Plots - Sanctuary
################################################################
pdf(file = "Figure3-D.pdf",   # The directory you want to save the file in
width = 9, # The width of the plot in inches
height = 3)
par(mar = c(1.5, 2, 2 , 1), oma = c(4, 4,1.5, 0))
par(mfrow = c(1, 5))
plot(x = data1c$time,
y = data1c$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-60, 60, 6),
xlab = "", ylab = "",
pch = 16, cex.axis=1.25)
segments(x0 = data1c$time, x1 = data1c$time,
y0 = data1c$lower_b, y1 = data1c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data2c$time,
y = data2c$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data2c$time, x1 = data2c$time,
y0 = data2c$lower_b, y1 = data2c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data3c$time,
y = data3c$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data3c$time, x1 = data3c$time,
y0 = data3c$lower_b, y1 = data3c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data4c$time,
y = data4c$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data4c$time, x1 = data4c$time,
y0 = data4c$lower_b, y1 = data4c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data5c$time,
y = data5c$estimate,
main = "",
ylim = c(-60, 60), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data5c$time, x1 = data5c$time,
y0 = data5c$lower_b, y1 = data5c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
mtext("Mahalanobis", line = -0.2, at = 0.11, outer = TRUE, cex = 1)
mtext("PS \n Matching", line = -1.8, at = 0.31, outer = TRUE, cex = 1)
mtext("PS \n Weighting", line = -1.8, at = 0.51, outer = TRUE, cex = 1)
mtext("CBPS \n Matching", line = -1.8, at = 0.71, outer = TRUE, cex = 1)
mtext("CBPS \n Weighting", line = -1.8, at = 0.91, outer = TRUE, cex = 1)
mtext("Estimated Effect of Sanctuary", side = 2, line = .3, at = 0.5,
outer = TRUE, cex = 1)
mtext(1, text = "Years relative to the administration of treatment", line = 1.2,
at = 0.5, outer = TRUE, cex = 1)
dev.off()
################################################################
# Robustness Checks 4: Civilians Killed by ELN (OMC)
################################################################
rm(list = ls())
load("RCh4-OMC_ResultsAlternativeModel_ELN.RData")
################################################################
# Preparing the data - Protests
################################################################
x1 <- summary(pe.results5.p2_maha5b)
x2 <- summary(pe.results5.p2_ps.m5)
x3 <- summary(pe.results5.p2_ps.w5)
x4 <- summary(pe.results5.p2_CBPSm5)
x5 <- summary(pe.results5.p2_CBPSw5)
data1a <- as.data.frame(x1$summary)
data1a <- data1a %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data1a <- data1a %>%
mutate(time = c(0, 1, 2, 3, 4))
data2a <- as.data.frame(x2$summary)
data2a <- data2a %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data2a <- data2a %>%
mutate(time = c(0, 1, 2, 3, 4))
data3a <- as.data.frame(x3$summary)
data3a <- data3a %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data3a <- data3a %>%
mutate(time = c(0, 1, 2, 3, 4))
data4a <- as.data.frame(x4$summary)
data4a <- data4a %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data4a <- data4a %>%
mutate(time = c(0, 1, 2, 3, 4))
data5a <- as.data.frame(x5$summary)
data5a <- data5a %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data5a <- data5a %>%
mutate(time = c(0, 1, 2, 3, 4))
################################################################
# Plots - Protest
################################################################
pdf(file = "Figure1-E.pdf",   # The directory you want to save the file in
width = 9, # The width of the plot in inches
height = 3)
par(mar = c(1.5, 2, 2 , 1), oma = c(4, 4,1.5, 0))
par(mfrow = c(1, 5))
plot(x = data1a$time,
y = data1a$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 20, 4),
xlab = "", ylab = "",
pch = 16, cex.axis=1.25)
segments(x0 = data1a$time, x1 = data1a$time,
y0 = data1a$lower_b, y1 = data1a$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data2a$time,
y = data2a$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 20, 4),
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data2a$time, x1 = data2a$time,
y0 = data2a$lower_b, y1 = data2a$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data3a$time,
y = data3a$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 20, 4),
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data3a$time, x1 = data3a$time,
y0 = data3a$lower_b, y1 = data3a$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data4a$time,
y = data4a$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 20, 4),
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data4a$time, x1 = data4a$time,
y0 = data4a$lower_b, y1 = data4a$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data5a$time,
y = data5a$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 20, 4),
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data5a$time, x1 = data5a$time,
y0 = data5a$lower_b, y1 = data5a$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
mtext("Mahalanobis", line = -0.2, at = 0.11, outer = TRUE, cex = 1)
mtext("PS \n Matching", line = -1.8, at = 0.31, outer = TRUE, cex = 1)
mtext("PS \n Weighting", line = -1.8, at = 0.51, outer = TRUE, cex = 1)
mtext("CBPS \n Matching", line = -1.8, at = 0.71, outer = TRUE, cex = 1)
mtext("CBPS \n Weighting", line = -1.8, at = 0.91, outer = TRUE, cex = 1)
mtext("Estimated Effect of Protest", side = 2, line = .3, at = 0.5,
outer = TRUE, cex = 1)
mtext(1, text = "Years relative to the administration of treatment", line = 1.2,
at = 0.5, outer = TRUE, cex = 1)
dev.off()
################################################################
# Preparing the data - Violent Self-Protection
################################################################
x1b <- summary(pe.results5.v2_maha5b)
x2b <- summary(pe.results5.v2_ps.m5)
x3b <- summary(pe.results5.v2_ps.w5)
x4b <- summary(pe.results5.v2_CBPSm5)
x5b <- summary(pe.results5.v2_CBPSw5)
data1b <- as.data.frame(x1b$summary)
data1b <- data1b %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data1b <- data1b %>%
mutate(time = c(0, 1, 2, 3, 4))
data2b <- as.data.frame(x2b$summary)
data2b <- data2b %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data2b <- data2b %>%
mutate(time = c(0, 1, 2, 3, 4))
data3b <- as.data.frame(x3b$summary)
data3b <- data3b %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data3b <- data3b %>%
mutate(time = c(0, 1, 2, 3, 4))
data4b <- as.data.frame(x4b$summary)
data4b <- data4b %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data4b <- data4b %>%
mutate(time = c(0, 1, 2, 3, 4))
data5b <- as.data.frame(x5b$summary)
data5b <- data5b %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data5b <- data5b %>%
mutate(time = c(0, 1, 2, 3, 4))
################################################################
# Plots - VSP
################################################################
pdf(file = "Figure2-E.pdf",   # The directory you want to save the file in
width = 9, # The width of the plot in inches
height = 3)
par(mar = c(1.5, 2, 2 , 1), oma = c(4, 4,1.5, 0))
par(mfrow = c(1, 5))
plot(x = data1b$time,
y = data1b$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 20, 4),
xlab = "", ylab = "",
pch = 16, cex.axis=1.25)
segments(x0 = data1b$time, x1 = data1b$time,
y0 = data1b$lower_b, y1 = data1b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data2b$time,
y = data2b$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data2b$time, x1 = data2b$time,
y0 = data2b$lower_b, y1 = data2b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data3b$time,
y = data3b$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data3b$time, x1 = data3b$time,
y0 = data3b$lower_b, y1 = data3b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data4b$time,
y = data4b$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data4b$time, x1 = data4b$time,
y0 = data4b$lower_b, y1 = data4b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data5b$time,
y = data5b$estimate,
main = "",
ylim = c(-20, 20), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data5b$time, x1 = data5b$time,
y0 = data5b$lower_b, y1 = data5b$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
mtext("Mahalanobis", line = -0.2, at = 0.11, outer = TRUE, cex = 1)
mtext("PS \n Matching", line = -1.8, at = 0.31, outer = TRUE, cex = 1)
mtext("PS \n Weighting", line = -1.8, at = 0.51, outer = TRUE, cex = 1)
mtext("CBPS \n Matching", line = -1.8, at = 0.71, outer = TRUE, cex = 1)
mtext("CBPS \n Weighting", line = -1.8, at = 0.91, outer = TRUE, cex = 1)
mtext("Estimated Effect of VSP", side = 2, line = .3, at = 0.5,
outer = TRUE, cex = 1)
mtext(1, text = "Years relative to the administration of treatment", line = 1.2,
at = 0.5, outer = TRUE, cex = 1)
dev.off()
################################################################
# Preparing the data - Sanctuary
################################################################
x1c <- summary(pe.results5.s2_maha5b)
x2c <- summary(pe.results5.s2_ps.m5)
x3c <- summary(pe.results5.s2_ps.w5)
x4c <- summary(pe.results5.s2_CBPSm5)
x5c <- summary(pe.results5.s2_CBPSw5)
data1c <- as.data.frame(x1c$summary)
data1c <- data1c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data1c <- data1c %>%
mutate(time = c(0, 1, 2, 3, 4))
data2c <- as.data.frame(x2c$summary)
data2c <- data2c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data2c <- data2c %>%
mutate(time = c(0, 1, 2, 3, 4))
data3c <- as.data.frame(x3c$summary)
data3c <- data3c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data3c <- data3c %>%
mutate(time = c(0, 1, 2, 3, 4))
data4c <- as.data.frame(x4c$summary)
data4c <- data4c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data4c <- data4c %>%
mutate(time = c(0, 1, 2, 3, 4))
data5c <- as.data.frame(x5c$summary)
data5c <- data5c %>%
rename(lower_b = "2.5%",
upper_b = "97.5%")
data5c <- data5c %>%
mutate(time = c(0, 1, 2, 3, 4))
################################################################
# Plots - Sanctuary
################################################################
pdf(file = "Figure3-E.pdf",   # The directory you want to save the file in
width = 9, # The width of the plot in inches
height = 3)
par(mar = c(1.5, 2, 2 , 1), oma = c(4, 4,1.5, 0))
par(mfrow = c(1, 5))
plot(x = data1c$time,
y = data1c$estimate,
main = "",
ylim = c(-20, 10), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), yaxp = c(-20, 10, 6),
xlab = "", ylab = "",
pch = 16, cex.axis=1.25)
segments(x0 = data1c$time, x1 = data1c$time,
y0 = data1c$lower_b, y1 = data1c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data2c$time,
y = data2c$estimate,
main = "",
ylim = c(-20, 10), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data2c$time, x1 = data2c$time,
y0 = data2c$lower_b, y1 = data2c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data3c$time,
y = data3c$estimate,
main = "",
ylim = c(-20, 10), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data3c$time, x1 = data3c$time,
y0 = data3c$lower_b, y1 = data3c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data4c$time,
y = data4c$estimate,
main = "",
ylim = c(-20, 10), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data4c$time, x1 = data4c$time,
y0 = data4c$lower_b, y1 = data4c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
plot(x = data5c$time,
y = data5c$estimate,
main = "",
ylim = c(-20, 10), xlim = c(-0.5, 4.5),
xaxp = c(0, 4, 4), xlab = "", ylab = "",
pch = 16, yaxt = "n", cex.axis=1.25)
segments(x0 = data5c$time, x1 = data5c$time,
y0 = data5c$lower_b, y1 = data5c$upper_b,
lwd = 1, col = "black")
abline(h = 0, lty = 2)
mtext("Mahalanobis", line = -0.2, at = 0.11, outer = TRUE, cex = 1)
mtext("PS \n Matching", line = -1.8, at = 0.31, outer = TRUE, cex = 1)
mtext("PS \n Weighting", line = -1.8, at = 0.51, outer = TRUE, cex = 1)
mtext("CBPS \n Matching", line = -1.8, at = 0.71, outer = TRUE, cex = 1)
mtext("CBPS \n Weighting", line = -1.8, at = 0.91, outer = TRUE, cex = 1)
mtext("Estimated Effect of Sanctuary", side = 2, line = .3, at = 0.5,
outer = TRUE, cex = 1)
mtext(1, text = "Years relative to the administration of treatment", line = 1.2,
at = 0.5, outer = TRUE, cex = 1)
dev.off()
rm(list = ls())
################################################################
# End of this part
################################################################
